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The supercritical series expansion of the survival probability for the one-dimensional contact pro- 
cess in heterogeneous and disordered lattices is used for the evaluation of the loci of critical points 
and critical exponents (5. The heterogeneity and disorder are modeled by considering binary regular 
and irregular lattices of nodes characterized by different recovery rates and identical transmission 
rates. Two analytical approaches based on Nested Pade approximants and Partial Differential ap- 
proximants were used in the case of expansions with respect to two variables (two recovery rates) for 
the evaluation of the critical values and critical exponents. The critical exponents in heterogeneous 
systems are very close to those for the homogeneous contact process thus confirming that the con- 
tact process in periodic heterogeneous environment belongs to the directed percolation universality 
class. The disordered systems, in contrast, seem to have continuously varying critical exponents. 
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I. INTRODUCTION 

In nonequilibrium statistical mechanics, phase tran- 
sitions have been identified and studied for some time 
now. Similar to equilibrium systems, these phase transi- 
tions fall into a number of different universality classes, 
one of which is the directed percolation (DP) universality 
class [Hll]- According to a conjecture by Janssen Q and 
Grassberger [3| , all absorbing state phase transitions with 
a scalar order parameter and no additional conservation 
laws are characterized by DP critical exponents. 

One of the questions that has not been answered con- 
clusively is how the introduction of quenched disorder 
affects the universal critical behavior of the DP class. Ac- 
cording to the Harris criterion Q , the critical exponents 
should change even for weak disorder. This has been in- 
vestigated in particular for the contact process (CP) Q 
which is one of the archetypical models of the DP univer- 
sality class. Currently, there are two alternative scenarios 
how the exponents change with introduction of disorder: 
according to results in Refs. 

0, i, S [M [HI , they change 
continuously with degree of disorder, while Vojta [l^ has 
suggested an abrupt change to the values in the strong 
disorder limit corresponding to the universality class of 
the random transverse Ising model. 

The above controversy can be addressed either by nu- 
merical or analytical methods. Among analytical meth- 
ods, time-dependent perturbation theory as introduced 
by Dickman and Jensen [13, 14] for homogeneous Id sys- 
tems gives the most accurate numerical values. Tech- 
nically, this has been done by using one-variable Pade 
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approximants in the analysis of the series for the survival 
probability. In this paper, we extend their approach to 
heterogeneous and disordered Id lattices and introduce 
the use of a Pade approximants similar to the Nested 
Pade approximants (NPA's) [H, [l^ in order to deal with 
two control parameters (two recovery rates), character- 
istic of binary lattices. A different approach based on 
partial differential approximants (PDA's) has been used 
by Dantas and Stilck in Ref. [13] , who applied the super- 
critical series expansion to study the crossover between 
the Id CP and the voter model thereby introducing 
a second control parameter to the perturbation theory. 
We also use PDA's in order to compare results from the 
two extrapolation methods. 



The aim of the paper is two-fold: first, we present 
and discuss the technical details of the supercritical se- 
ries expansions in the case of two variables, i.e. of two 
different recovery rates characteristic of nodes of two dif- 
ferent types, which are different from the well-studied 
one-parameter case [H, [l3]- Second, we investigate the 
range of applicability and the effects of variations of our 
procedure based on NPA's on the estimates of critical 
values and exponents and compare these results to those 
obtained by employing PDA's. 



The structure of the paper is as follows: we introduce 
the CP in Sec. |TT1 The supercritical series expansion, the 
analysis of the resulting two- variable series and the con- 
figurational averaging of the order parameter for disor- 
dered environments arc discussed in Sec. IIIII In Sec. IIVI 
we present the results in terms of phase diagrams and 
the critical exponent /S for the different systems, which 
we then discuss in the last Sec. IVl 
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II. MODEL 

The CP, originally introduced by Harris [3], is a spa- 
tial SIS (Susceptible-Infected-Susceptible)-model for the 
spread of epidemics through a network. In a network, 
usually taken to be a hyper-cubic lattice, the nodes or 
sites can be in one of two states: "infected" or "suscep- 
tible". A susceptible/ vacant site i can be infected by a 
neighboring infected/occupied site j, while an occupied 
site k can spontaneously recover and become susceptible 
again. These processes occur with rates Xij/z and fik, 
respectively, where z is the number of nearest neighbors 
in the lattice. 

In all dimensions, the CP undergoes a nonequilibrium 
phase transition into an absorbing state which does not 
allow any further time evolution. For a fixed set of trans- 
mission rates, this occurs when the recovery rates become 
smaller than a set of critical values. At this point, the 
survival probability, that is the probability that process 
will not enter the absorbing state (in the thermodynamic 
limit), becomes greater than zero. Several observables 
can be identified that describe the critical state of the 
CP, such as the average density of infected sites, p (t), or 
the survival probability up to time t after starting from a 
single seed, Pg {t)- As t — *■ oo, in a homogeneous system 
(Xij — 1 and fik ~ A*), close to the critical point these 
quantities are expected to behave like ^ 

lim p{t) = pstat ~ A'^ lim Psit) = Poo^ A^. (1) 

t — >QO t — >QO 

Here P is the critical exponent associated with stationary 
behavior of the order parameter and A = — M > 0, 
where fie is the critical point for the control parameter 
recovery rate /x. 

Each site i of the system of Ng nodes obeying the rules 
of the CP can be in two states: \ai) where ai = or 
CTi = 1 so that a microstate of the system with Ns sites 
can be written as 

|a) = |ai)0...®|a^J=(g)K) = (g)f ^^'^^ j 

This representation ensures that the microstate vectors 
form a 2^= -dimensional orthonormal basis. The state of 
the system at time t is 

|P(i))=^P(a,i)|a) (2) 

where P{a,t) = {a\P{t)) is the probability that the sys- 
tem is found in microstate |cr). The time-evolution of the 
state of the system is governed by the master equation 

d,\P{t)) ^ L\P{t)) (3) 

where L is the Liouville operator whose non-zero off- 
diagonal elements in this basis are the transition rates 
between microstates that for the CP differ in their oc- 
cupation number by one. This operator describes the 



probability flow between different microstates and is thus 
represented by a stochastic matrix in which all the diago- 
nal elements are the sums of the off-diagonal elements in 
the corresponding columns taken with the opposite sign. 

In a formalism introduced by Doi [19] and Peliti [2^ 
for stochastic systems, "annihilation" and "creation" op- 
erators on site i, and aj, respectively, are defined, such 
that |(Ti) = (Ti \ai - 1) and aJ |(Tj) = (1 - oi) \ai -f 1) 
(e.g. a-|0) = 0) and that they obey hard-core bosonic 
commutation relations. 

For simplicity, we assume all the transmission rates to 
be the same and define the time scale by setting Ay = 1 . 
The recovery rates for binary systems are characterized 
by one of two values pi = {\xa-,Pb\- For such models, 
the operator L in the onc-dimensional case reads as 

L = pW^V, (4) 
where the operators p W and V are 

= ^^/fl(l-at)a, , (5) 

i 

The parameter fi is introduced for convenience and 
~ 1 with both recovery rates being close to the 
homogeneous critical point and thus Pi 1- The op- 
erator V creates offspring at the nearest neighbors of an 
occupied sites, while W destroys occupation at sites. The 
above formalism is useful for the supercritical series ex- 
pansion in fi described below. 

III. ANALYSIS 

A. Supercritical Series Expansion 

The supercritical series expansion is a perturbation se- 
ries for the ultimate survival probability Poo which is 
taken to be the order parameter, with Poo > being 
characteristic of the CP in the active phase. To probe 
the long-time limit of the system, the Laplace transform 
of the probability state vector is taken, 

poo 

|P(.))= / dte~^'\P{t)) , (7) 
Jo 

so that a standard identity of Laplace-transform theory, 
limt^oo fit) — lims^o s /(s), can be used. The ultimate 
survival probability is then given by 

Poo = ^rn(l-s{0\P{s))) . (8) 

Inserting the formal solution |P(t)) e"^*|P(0)) of 
Eq. dni) into Eq. ^ resuhs in 

\Pis))^{s~pW-V)-'\P{0)) . (9) 
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We can then formally expand the operator on the right 
hand-side of Eq. ([H]) in a Taylor series with respect to 
the small parameter fj, thus obtaining the following su- 
percritical expansion 

|P(.)) = \Pois))+^,\P,is))+f,'\P,is)) + ...{10) 

where the vectors |P„(s)) obey the following recursion 
relations: 

\Po{s)) = {s-V)-'\PiO)) (11) 
\Pn{s)) - (s - Vr'W\Pn-iis)) for n > 1 . (12) 

The action of the operator on a given configuration 
can be straightforwardly computed using its definition 
given by Eq. (O, i.e. 

m 



The summation is taken over all m occupied sites in state 
a and = ajo-), i-e. the ct* and a have the same 
occupation except for site i being occupied in state a 
and vacant in a^. 

The action of the operator (s—V)^^ on a given configu- 
ration, \<t), can be computed with the use of the following 
identity, 



{s - V)-^ = s-^ + s~\s - V)-^V , (14) 



so that 



{s-V)-'\<7) = 



1 



s + qi/2 + 52 



\(j) + {s-vy 



-, 91 92 



(15) 



where the sums represent the action of the operator V 
on the state a with vacant sites of two types: sites which 
have one (first sum in Eq. (|15p is taken over a number 
qi of such sites) or two (second sum in Eq. p5|) is taken 
over of such sites) occupied nearest neighbors. The 
vectors and represent the states in which the 

formerly vacant sites i of the first and second type, re- 
spectively, are now occupied. To go further we should use 
the recursive nature of the operator (s — V)~^ and sub- 
stitute its representation given by Eq. into Eq. 
Such a procedure can generate an infinite number of new 
configurations. However, when we calculate the survival 
probability perturbatively up to a given order N in fi 
for the initial condition of a single occupied site, it is 
only necessary to retain states with up to N occupied 
sites. This is due the fact that the annihilation opera- 
tor W in an expansion up to order /i^ acts N times on 
any generated state, after which remaining states will be 
projected onto the absorbing state, thereby projecting 
out any states with more than N occupied sites. Follow- 
ing this procedure, we can perturbatively calculate the 
survival probability Poo and thus find the critical point 
where the survival probability becomes zero. In order 
to obtain good numerical estimates of this critical value 
and compute critical exponents, it is necessary to employ 
numerical methods such as Fade approximants plj . 



B. Nested Pade Approximants: Critical Values and 
Exponent 

For systems with two different recovery rates, iia and 
/is, the survival probability (see Eq. ([8])) expanded in 
series is a polynomial in these two variables, 

N n 

Poo (ma, /is) = 1 - E E ^a"" ■ (16) 

n— 1 m— 

The critical line /i^-* — /i^'(/i^') that separates the ab- 
sorbing state from the active state is a solution to the 
equation Poo(/iA, /^s) — corresponding to the small- 
est (real) root. In practice, just finding the roots of the 
polynomial - the truncation of an infinite series at fi- 
nite order - does not produce very good estimates of 
the critical values. Better results are obtained by us- 
ing d-log Fade approximants (2^ : given an expansion in 
one variable. Poo (m) , up to order TV, the Fade approxi- 
mant [2l|, [23| of the a series for 9^ In P^o (/i) is formed. 
Technically, this is done by expanding the denominator 
of 9^ In Poo (/i) = d^P{n) / P{ii) up to order TV — 1 and 
thus obtaining a polynomial /Ar-i(/i) for this fraction, 
the Fade approximant of which is then constructed. The 
first positive (real) pole and its residue then provide good 
estimates of the critical value and the critical exponent 
of Poo, respectively. The extension of this approach to 
two variables, however, is not straightforward - there is a 
number of different multi- variable generalizations (2^ . [25j 
of the one- variable Fade approximation. In this work, we 
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employ a scheme similar to the NPA's [H, [T^ in which 
we in turn form one- variable Pade approximants with re- 
spect to the two variables. 

To this end, we transform the variables (/ia, ^J■B) to the 
following three more convenient sets, Tl, T2 and T3. The 
first transformation, Tl, is symmetric so that the values 
{flA,^J■B) are replaced by {ft, 5) where ft = {^a + 
and S = {fiA — fJ'B)/'2- The transformations T2 and T3 
are asymmetric with (^a, ^J■B) replaced by either {p, = 
fiA,S ^ ^lA-^J^B) or {fL = fiB.S = ^iB-fJ-A), respectively. 
Expanding ^ InPaoift, S) = dp_Poo/ Pooift, S) up to order 
N — 1 in ft and S, we obtain 



9/i 



(17) 



where N— 1 — n = J + K. As always with Fade ap- 
proximants, we have the freedom to choose J and K for 
a given n. For even N — 1 — n, we use diagonal Fade 
approximants [J, J] with J — K — {N — 1 — n)/2, while 
for odd iV — 1 — n, we form approximants [J — 1, J] with 
J= {N -n)/2. 

Then, in turn, we form the Fade approximant with 
respect to /2, 



i9m In ^o< 



ELornjS) fl- 

I + En=l rL+n{6) ft- 



+ 0(/2^) (19) 



We then form the Fade approximants of the coefficients. 



fn{S) 



Em=0 9™'^" 



1 + Em=l qj+mS"" 



(18) 



where iV — 1 = L + Af . A more graphical representation 
of the scheme is the following: 



{ao.o + aa^iS + . . . + aQ^N-i6 ) + {ai,o + . . . ai,N-i S ) + . . .ajq-ifl 
^ ^ ' ' 

[[{N ~ 1)/2J, \{N - l)/2]], + [[{N - 3)/2j, \{N - l)/2]], /I + . . . a^_i,oM 



[L(A^-l)/2j,r(7V-l)/2]]_+0(/2^) 
I 



+ 0(/2^) (20) 
^-i+0(a2^) (21) 
(22) 



Here we have denoted the formation of the Fade approxi- 
mant with respect to a variable x with numerator degree 
N and denominator degree M as [TV, M] ^ ( [-J and [•] 
are the floor and ceiling functions, respectively). 

Thus for any given 5, we can find the corresponding 
pole of djj. In Poo (m, <5), which is then taken to as the crit- 
ical value ft^'^\6), yielding a point on the critical line, 

{^^^^ = /iw,^!,^) = m(^) + ^), {^^^^ = m'^' = M^^)) 

or (/^^^ — ft^"^^ + 5,^^^ = ft^'^^ — 5), depending on the 
initial transformation. 

It turns out that occasionally the first positive real 
roots are unphysical ones that appear before the phys- 
ical solution. However, these roots are very closely 
matched by roots of the numerator of the Fade approx- 
imant of 9p lnPoo(MAj (5), so that these two cancel each 
other, leaving the physical root as the solution. In or- 
der to extract unphysical roots a further parameter 7 
has been introduced. Two roots, Xi and X2, of the nu- 
merator, n{x), and the denominator, d{x), respectively, 
i.e. {x — x{)n{x)/ {{x — X2)d{x)), are considered to be the 
same value and cancel if |a:i — X2\ < 7. All results pre- 
sented below are obtained by setting 7 = 10""^. 

In order to evaluate the stability of a certain pole, sev- 
eral approximants are formed with respect to ft close to 
the diagonal approximant, e.g. for even A'^ — 1 = 2K, we 



compute [K, K],[K - 1, K],[K, K - 1], . . .[K - 2, K - 2] 
and take the average over these poles. Once we obtained 
the critical value, the critical exponent (3 associated with 
the order parameter Poo can be found as well, as it is just 
the residue at the pole fi^'^\ Again, the average over the 
residues for different Fade approximants is taken. As the 
error that we could extract from employing this method 
is small and does not take into account the inherent er- 
ror in the series expansion, we perform this averaging to 
minimize the effects due to a particular choice of approx- 
imant and use the standard deviation only to evaluate 
the numerical stability of a pole. 



C. Partial Differential Approximants: Critical 
Values 

Another method for estimating critical values given a 
finite two-variable series are the FDA's originally devel- 
oped by Fisher [2^ in order to investigate multicritical 
points. 

The starting point of this method is that one 
is given a finite two-variable polynomial, F{x,y) = 
E(j: j)es /y '^*^"'' that approximates a true function 
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f{x,y), that is expected to have the following form, 

fix, y) = {{x - X,) + aiy ~ yc)f<f>ix, y) (23) 

where 4){x, y) is some general function with 0(a;c, yc) 7^ 0. 
The set S is a so-called label set that contains the 
pairs of powers of x and y only if ^ 0, i.e. 

G S if F{x,y) has a non-zero term of order x^y^ . 
For fixed label sets L, M and N, it is possible to find 
polynomials Pi.{x,y) = J2(z,j)ei,Pijx'y' ^ Qm{x,y) = 

T.{;,j)e_m1^i^'y^ a-iid ^N(a;,y) = E(»j)eN 2;*^^ such 
that they satisfy the defining equation 

dF dF 
Pi.{x,y)F = Q^{x,y)-^ RN{x,y)-^ + Ei<i{x,yl2A) 

with E'k = '^/(i j)^vi^ijX^y^ denoting a sum of non- 
zero terms whose powers are not in the matching set K. 
This matching set defines for which powers x'^y^ Eq. ([M]) 
should hold exactly with ~ 0, while for {m,n) ^ K, 
the values of Cmn are allowed to be non-zero. In order 
for Eq. to be a solvable set of linear equations, the 
label sets must obey the constraint that the label sets L, 
M and N must together contain one more element than 
the matching set K because of the conventional choice of 
Poo = 1- 

Once the polynomials Pl, Qm and Pn are found, e.g. 
by using an algorithm proposed by Styer [23| , they can be 
used to find an estimate for the line of critical points by 
the method of characteristics (e.g. see [111)- According to 
this method, consider a single curve of points which only 
depends on a single parameter r, x(t) — {x(T),y(T)). 
The rate of change of F{x{T),y{T)) along this line is 
dF/dr = {dF/dx)idx/dT) + {dF/dy){dy/dT). The sur- 
vival probability Poo , for which we are going to apply this 
method, is zero along the critical line. Thus, considering 
the case where F{x{T),y{T)) = 0, it can be seen that the 
curve described by the equations 

^ - Qm{x{r),y{r)) (25) 

g = PN(x(r),y(r)) (26) 

and substituted in Eq. ([24|) yields the relation = 
{dF/dx){dx/dT) + {dF/dy){dy/dT). Together with a 
suitable initial condition, this curve is therefore equiv- 
alent to the critical line as F does not change along the 
curve x(r). This suitable initial condition has to be a 
known point on the critical line: in our case, this is 
the critical point of the homogeneous system at which 
X = y ^ Xc- 

D. Configurational Averaging 

The schemes above are straightforwardly applied to 
heterogeneous topologically ordered systems. In topo- 
logically disordered systems, the survival probability has 



to be averaged over different realizations of disorder. For 
concreteness, we consider disorder only in recovery rates 
fii on different nodes i assuming that are independent 
random variables distributed according to the probabil- 
ity distribution functions p{fJ-i). A configurationally av- 
eraged survival probability is then given by the following 
expression, (Poo) = / Poo(Aii) Di P (Mi)d/ij- For simplic- 
ity, we consider a bimodal distribution of recovery rates, 

p(m») =P'5(^i - ha) + (1 ~p)S{fii ~ ^ib) , (27) 

i.e. the nodes A (hosts), characterized by recovery rate 
fiA, of concentration p and B (impurities) of concentra- 
tion 1 ~p are randomly and independently placed on the 
sites of a regular chain. 

Using the series expansion of order N for Poo is equiv- 
alent to considering the CP on the finite chain of length 
2A^ — 1, i.e. for a given value oi n < N all states on 
2n — 1 sites with the origin at its center with at most n 
sites occupied contribute, so that 

N n 

(Poo) = 1 + E E ^'B A^r™ (28) 

n— 1 m— 

where 

2n-l 

(c„™)=E E (l-p)V""'~'c„„(c) (29) 

k=0 {c:NB{c)=k} 

with the second sum running over all disorder configura- 
tions c — (/i_„+i, . . . , /io, ■ • ■ ) A*n-i) that have a certain 
number Nb{c) of impurity sites B. The values of Cnm 
are the coefficients in the expansion of the survival prob- 
ability for a particular disorder configuration. The factor 
(1 — p)'^p2n-i-A; g^gj^jg fj-Qm the fact that the probabil- 
ity of having a particular disorder configuration is just 
the product of the probabilities of any site being either 
fiA or fiB , drawn from the bimodal distribution given by 
Eq. 1271). 

The memory requirements in calculations of the co- 
efficients in the series expansions impose a restriction, 
N < A'max, on the highest order of expansion for disor- 
dered lattices, iVmax = 19, which is rather lower than for 
homogeneous [l^ and heterogeneous cases, A^max = 24. 
The exact configurational averaging discussed above, ex- 
ploiting symmetry about the origin and under the ex- 
change of iJ,A and /is , is not possible for such high orders 
due to computational cost when dealing with a very large 
number (of 0(2^^"""')) of configurations. We have been 
able to undertake the exact configurational averaging up 
to order A^cmax = 12. 

For higher orders, A^cmax < A^ < A^max, a config- 
urationally averaged survival probability is calculated 
approximately by only including disorder configurations 
that have no more than a certain number of impurities 
in the averaging. Assuming that each coefficient is of 
the order of the coefficient in a homogeneous system, 
Cnm(c) ^ c„, then each coefficient (c„m) will remain of 



the same order of magnitude if we choose a maximum 
number of impurities k^a,-^{n,p) in 

k 

{Cnm) - Cn E " ^) (1 - P) V"-^-'= (30) 

such that the sum is close to unity. For (1 — p) <C 1 
and large n ^ 1, one can choose fcmax ^ 2n — 1, e.g. 
for fc^a,,(19,0.96) = 2, ^^^0" ('T')(l " V"~^-'= = 
0.817. For lower orders, the weight of the configura- 
tions that are dropped decreases for constant fcmax, but, 
because it becomes computationally feasible, we choose 
to increase fcmax by one for each lower order by letting 
fcmax(-^max — hp) =1 + 2. This way, the lower the order 
the closer the approximate configurational averaging is 
to the exact one. With fcinax(^max — hp) = « + 2, the 
averaging becomes exact from order n < 7V/3 + 1. 

IV. RESULTS 

In this section, we will present the results for criti- 
cal values and critical exponents obtained by the anal- 
ysis described above and applied to different systems. 
These results come in the form of phase diagrams in 
which the critical points are plotted in the rate-space 
plane {^a, Mb) or as plots of the critical exponent /? as a 
function of the critical rate ^J^\ ■ 

In Ref. 29], it has been demonstrated that for the Id 
CP the line of critical points close to the homogeneous 
critical point, (//a = McMs = Mc) is well described by 
the relation 

^exp(E[ln(/.,)]) , (31) 

where E[-] denotes the expectation value with respect to 
the distribution of the recovery rates /i^. In this work, 
we will compare the critical points to the line given by 
Eq. pip to examine how far from the homogeneous crit- 
ical point the relationship describes the critical line well. 

In all figures below, we only keep critical points and 
critical exponents that have a standard deviations from 
the mean values of less than 0.001 and 0.005, respectively, 
after averaging over the Fade approximants as described 
in Sec. HITB] 

A. Periodic Lattices 

First, we analyzed periodic lattices, i.e. the CP in sys- 
tems with a repeating pattern of nodes characterized by 
the recovery rates and ^b- Thinking of a periodic 
system in terms of unit cells that are repeated through- 
out the lattice and denoting a site in this unit cell that 
has recovery rate y,A {f^s) as an A{B)-site, the three Id- 
lattices AB, AAB and AABB are considered. As can 
be easily seen from Eq. ([?!]) ■ the critical lines of AB and 
AABB should coincide, at least sufficiently close to the 
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FIG. 1: (Color online) Comparison of critical points (a) and 
critical exponents (b) of the AB system obtained from series 
expansion for different orders A'^: = 10 (o), = 17 (o) and 
N = 2A (-I-). The symmetric transformation Tl was used. 

The dashed lines ( ) are given by Eq. (|3H) in (a) and by 

P at the homogeneous critical point, /3 = 0.2769, from series 
expansions 14] for the CP, in (b). 
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FIG. 2: (Color online) Comparison of the critical points (a) 
and critical exponents (b) of the AB system obtained from 
series expansions up to order A'' = 24 with different transfor- 
mations, Tl (o) and T2 (-I-). The dashed lines ( ) are given 

by Eq. ([31} in (a) and by /3 at the homogeneous critical point, 
/? = 0.2769, from series expansions TS] for the CP, in (b). 



homogeneous critical point. Therefore, we will only con- 
sider AB and AAB in detail, except for a comparison of 
the stability of the critical values away from the homo- 
geneous critical point for AB and AABB, which will be 
the subject of the last part of this section. 
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FIG. 4: (Color online) Comparison of critical values (a) and 
critical exponents (b) of the AAB lattice obtained from series 
expansions up to overall order N — 24 which are cut off at 
order A/ — 12 (+) or in which all the terms up to order 
M = 23 (o) are retained. Transformation T2 was used to 

obtain these critical points. The dashed lines ( ) are given 

by Eq. (|3ip in (a) and by f3 at the homogeneous critical point, 
/3 — 0.2769, from series expansions 14] for the CP, in (b). 



FIG. 5; (Color online) Comparison of critical values (a) and 
critical exponents (b) of the AB lattice obtained from series 
expansions up to overall order A*' = 24 which are cut off at 
order M = 12 (+) or in which all the terms up to order M = 
23 (o) are retained. Transformations T2 and Tl, respectively, 
were used to obtain these critical points. The dashed lines 

( ) are given by Eq. ((31} in (a) and by /3 at the homogeneous 

critical point, /3 = 0.2769, from series expansions fT3] for the 
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FIG. 6: (Color online) Comparison of phase diagram for the 
AB system obtained from series expansion up to order A*' = 12 
(o) and for the AABB chain obtained from series expansion 

up to order A'' = 24 (+). The dashed line ( ) is given by 

Eq. ((311) • 



1. Results obtained by NPA's 

The series expansion for all three systems has been cal- 
culated up to order N — 24. In order to evaluate how the 
estimates provided by the analysis described in Sec. lIIIBl 
behave with the order of the expansion, in Figs. [T] (a) 
and (b) the critical line and the critical exponent for 
the AB-lattice are shown for four different values of N, 



iV = {10, 17, 24}. The points in Figs. (D (a) and (b) were 
obtained after using transformation Tl. Most notably, 
as seen in Fig. (TJ the range of reliable critical points in- 
creases with increasing order N. Of course, the accuracy 
of the prediction also increases with increasing N. How- 
ever, in general it can be said that even for low orders of 
the expansion the critical values are of good accuracy. 
The critical exponents are more sensitive to the order 
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FIG. 7: (Color online) Critical lines in log-log scale for the 
AB lattice (o) and the AAB lattice (+ and <>) as obtained 

by PDA's: The dashed lines ( ) and (•• — ••) are given by 

Eq. (|5T|| for the AB- and the AAB-sjsiera^ respectively. 



of expansion - their behavior with N is shown in the 
Fig. [1] (b). Not surprisingly, we find that with increas- 
ing order the critical exponents come closer to the best 
known value for the homogeneous CP from series expan- 
sions, /3 ~ 0.2769 [3| - at least close to the homogeneous 
critical point. Further away from this point, all the crit- 
ical exponents for all orders fluctuate between (3 — 0.275 
and (3 = 0.28, the range in ha of reliable exponents co- 
inciding with the range for the critical points. 

In an attempt to extend the range of applicability of 
the series expansion to locate critical points, the linear 
transformations described in Sec. IIII Bl are applied to the 
expansion variables. These transformations change the 
magnitude of the new variables /2 and 5 and are therefore 
expected to provide extensions in different regions of the 
phase diagram. First, the transformations Tl and T2 
are in employed in the AB system (T3 being related by 
symmetry to T2). Fig. [2] shows that for this particular 
lattice the transformation Tl performs far better than 
T2 in giving stable critical points and critical exponents. 
The results of the transformations for the asymmetric 
system AAB are compared in Figs. [3] (a) and (b). Here, 
two things are to note: transformation T2 extends the 
critical line further in the < Mc direction than Tl and 
T3 does in general poorly compared to the other two. 

We also investigate the effects of different orders 
of Fade approximants for p. and S in the analy- 
sis as described in Sec. IIII Bl In order to do this, 
the steps given by Eqs. ((20)) and pT|) are followed, 
but then, instead of taking the Fade approximant 
[[(TV- 1)/2J, [(iV- l)/2]]- as shown in Eq. the 
series in Eq. (|2T|) is truncated at order AI < N — 1 
in fl. This means that the coefficients of the terms /l" 
with n < M that remain in the series are Fade approx- 
imants in 5 formed from polynomials of degree larger 



than A'' — 1 — M. Therefore, the extrapolation provided 
by the Fade approximants of these coefficients can be ex- 
pected to be more accurate than for polynomials of lower 
degree, which is the case for the coefficients of /I" for 
n > M. Then, the Fade approximant [[M/2\, [M/2]]^ 
is taken. The results of such a cut-off is that the range 
of the convergence of the series is improved allowing for 
well behaved poles further away from the homogeneous 
critical point. The extension of the range for the phase- 
separation line achieved by this procedure can be seen in 
Figs, m (a) andm (b) for the AAB-system. For this sys- 
tem, the transformation T2 significantly extends the line 
of critical points into the region /z)^ < /Ltc and /z^ > ^c- 
In Fig. m (b), it can be seen that the critical exponents 
deviate from the value /3 = 0.2769 by less than 1.2%. 

The fact that this procedure of truncating in orders of 
expansion in /l performs well in this particular case can 
be understood in the following manner. For the transfor- 
mation T2, p, = fiA and S = fiB — ft = Mb — Ma, leading 
to small values of the variable (p,) in which the expansion 
is cut, while S is very large for ^a < Mc and fiB > Mc- 
Therefore, we can expect the effects of dropping terms 
in p to be small while the benefit of only using Fade 
approximants in 6 formed from high-order polynomials 
to be large. This explains the sizeable extension of the 
region in which the poles of the series are well-behaved. 

In contrast to the AAB system, the procedure de- 
scribed above does not produce much extension of the 
critical line for the AB lattice, as can be seen Fig. [5] 
there are only a few points obtained by using transforma- 
tion T2 further away from the homogeneous point than 
the critical points by the regular analysis with Tl. 

It is worth noting that the deviation of critical line 
from the line given by Eq. (|3ip appreciably increases 
away from the homogeneous critical point (see Fig. [5]). 
It shows that Eq. (j3ip is only a lowest-order approxima- 
tion in relation to the homogeneous critical point. 

For the AABB system, it is found that, in general, 
there is a smaller range of the critical values and crit- 
ical exponents compared to the situation for AB and 
AAB chains. In order to investigate this difference in 
range of stability of the critical points and exponents, 
the AABB critical line is compared with that for the 
AB lattice, which should follow the same curve accord- 
ing to Eq. (|3ip . It is found that the difference of the 
stability of the critical values seems to be a result of the 
number of unit cells that the series expansion takes into 
account. For order N = 24, six unit cells on either side 
of the origin affect the series expansion for AABB while 
the number is twice that for AB. If order iV = 12 ex- 
pansion for the ^i?-lattice, which also only takes six unit 
cells into account, is compared with the order N — 
expansion for the AABB chain, then we find that their 
ranges are very similar. This is shown in Fig. [6] where the 
critical points start to fluctuate wildly for both system 
at the same point in the phase diagrams (not presented 
in Fig.©. 



2. Results obtained by PDA 's 



As described in Sec. IIII CI by applying the PDA 
method to the PooifJ'Aj ^J■A) we can also compute a line 
of critical points. The results obtained by this method 
depend on the choices of the label sets, L, M, N for 
the polynomials Pi^{x, y), Qm{x, y) and Rjssix, y) and the 
matching set K. For all our analyses employing PDA's, 
we use either label sets M, N and K that have trian- 
gular or rectangular form. The label set L is then cho- 
sen to mimic the form of the others while at the same 
time making sure that the number of elements in all four 
sets satisfies the constraint that is imposed on them (see 
Sec. IIII Cp . By "triangular" label set M we mean that 
e M if < (i + j) < M^ax while "rectangular" 
refers to a set M for which (i, j) G M if < i < Mmax,i 
and < j < M„iax.j with some integers M„iax-, Mmax,i 
and Mmax,j- From here on, we will refer to a particular 
choice of the label sets L, M, N and matching set K as 
an input set. 

Often, different input sets produce critical lines of vary- 
ing extent in (/iyi, /iB)-space. In Fig. [71 we show the re- 
sults from one input set for the AB lattice and from two 
input sets for the AAB lattice in comparison to the lines 
given by Eq. (pij) . The input set for the AB lattice are 
triangular while the ones for the AAB system are rect- 
angular and triangular. The critical lines obtained by 
PDA's are presented in log-log scale because they extend 
to a rather wider range than those calculated by NPA's 
(cf. Figs. [U [3] and [7]). One can clearly see that around 
the homogeneous critical point the lines from the series 
and described by Eq. (|3ip agree very well while devia- 
tions develop further away. It can be seen in Fig. [7] that 
the critical lines for triangular and rectangular input sets 
in the case of AAB chain coincide in the vicinity of the 
critical point with a consistent tendency to be above the 
prediction given by Eq. ([?T|) . There is a point though 
where these two approximate curves for the critical line 
diverge: from there, the estimates given by the PDA's 
can no longer be considered reliable (cf. the behavior of 
the curves marked by -I- and o in Fig. [7] for fiA ^1). 



B. Disordered Contact Process with Bimodal 
Distribution 

1. Results obtained by NPA's 

For the disordered CP, in which the recovery rates are 
drawn from the bimodal distribution given by Eq. (j27|) . 
the same analysis as for the periodic lattices is carried 
out, with the only difference that the survival probabil- 
ity is configurationally averaged. Up to order N — 12, 
the expansion of the survival probability can be configu- 
rationally averaged exactly for any value of p. The results 
for such systems characterized by p = 0.5 and p — 0.7 
are shown in Figs. [5] and [HI For both values of p, the 
line of critical points and the critical exponents for the 
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FIG. 8: (Color online) Comparison of critical values (a) and 
critical exponents (b) of the disordered system with p — 0.5 
obtained from series expansions up to order A'^ = 12 using 

transformations Tl (o) and T2 (-I-). The dashed lines ( ) 

denote (a) the curve given by Eq. pip and (b) the value for 
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FIG. 9: (Color online) Comparison of critical values (a) and 
critical exponents (b) of the disordered system with p = 0.7 
obtained from series expansions up to order = 12: Tl 
(o), T2 (-I-) and T3 (O). The dashed lines (---) denote (a) 
the curve given by Eq. (|3ip and (b) the value for the critical 
exponent at the homogeneous critical point for the CP from 
series expansion, /3 — 0.2769 [T^ . 



the two different linear transformations, Tl and T2, are 
compared. 

Similar to the heterogeneous AB lattice, we find that 
for the disordered system with p — 0.5, transforma- 
tion Tl extends the line of critical points furthest (see 
Fig. [5]) . The critical points are again well described by 
Eq. (|3ip . The critical exponents, however, behave differ- 
ently for the disordered case than for the heterogeneous 




FIG. 10: (Color online) Comparison of critical values (a) 
and critical exponents (b) of the disordered Id lattice with 
p = 0.96 obtained from series expansions up to order = 12 
using exact configurational averaging (o) and approximate 
configurational averaging (+) with fcmax(12,p) = 2. In both 
expansions Tl was used. The dashed lines ( ) are given by 
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FIG. 11: (Color online) Comparison of critical values (a) 
and critical exponents (b) of the disordered Id lattice with 
p = 0.96 obtained from series expansions of the approxi- 
mately configurationally averaged survival probability up to 
order N = IQ using the transformations Tl (o) and T2 

(+). The dashed lines ( ) denote (a) the curve given by 

Eq. (|31|l and (b) the value for the critical exponent at the 
homogeneous critical point for the CP from series expansion, 
/3 = 0.2769 JJ]. 



lattices: while in the latter case they fluctuate up imme- 
diately away from homogeneous point ix^ but then fluc- 
tuate down again, in the former system, they almost lin- 
early increase away from the homogeneous critical point. 

Fig. [H] shows the results for the disordered lattice with 
p = 0.7. This system behaves more like the heteroge- 



FIG. 12: (Color online) Critical lines in log-linear scale for the 
disordered lattice with p = 0.5 (o), p — 0.7 (+) and p — 0.96 

(o) as obtained by PDA's: The three lines ( — ), ( ) and 

(■■ — ■■) are given by Eq. (|3H) for the three disordered systems 
characterized by the values p = 0.5, p — 0.7 and p = 0.96, 
respectively. 



neous AAB lattice with transformation T2 extending the 
critical line furthest to the left, for iia < /^c and Tl fur- 
thest to the right, ha > Mc- As usual, the points obtained 
from T3 lie in the middle and thus do not extend the 
critical line in any direction. The critical exponent (3 dis- 
plays very similar behavior as for p = 0.5, mostly mono- 
tonically increasing away from the homogeneous critical 
point. 

As noted in Sec. IIIIDi the exact configurational aver- 
aging of the survival probability, (Poo), becomes rapidly 
computationally very demanding with order of expansion 
and is only possible for N < 12. Therefore, the approx- 
imate scheme for averaging has to be applied for orders 
N > 12. Such a scheme has been described in Sec. IIII P] 
In order to test the reliability of this scheme, we ana- 
lyze the disordered systems characterized by p = 0.96, 
for which both exact and approximate averaging, with 
^max(12,p) = 2, are possible. The comparison between 
the two configurational averages are shown in Fig. [TO] (a) 
and (b). It can be seen that the critical points agree 
very well for the two configurational averages where they 
both produce stable critical points, showing that the ap- 
proximation scheme is indeed faithful for small disorder 
concentrations. 

For a system with small concentration of B sites, 
e.g. for p > 0.96, the approximate scheme described in 
Sec. IIII Pj can be applied for configurational averaging up 
to order TV = 19. Below, we present results for p = 0.96 
with fci„ax(19,p) = 2. Figs.[Tl](a) and (b) show the crit- 
ical values and critical exponents for this case. Clearly, 
the transformation T2 works best for this system as is to 
be expected from its performance for the heterogeneous 
periodic AAB lattice. The critical points from transfor- 
mation T3 cover less range than Tl or T2, so we left 
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these points out of Fig. [TTJ Deviation of the critical Hne 
from the curve given by Eq. (|3ip can be seen for < Mc- 
The critical exponent (3 increases monotonically with in- 
creasing value of away from the homogeneous critical 
point, pi A > /ic) while for ^a < significant fluctua- 
tions can be seen. 

It should be mentioned that the critical line is only 
marginally extended by the procedure of cutting the se- 
ries in different orders for fl and S described in the pre- 
vious section. 



2. Results obtained PDA 's 

Using similar input sets as for the heterogeneous sys- 
tems (as described in Sec. lIV A"2t . we apply PDA's to the 
disordered systems as well. Generally, smaller input sets, 
with fewer elements in the label sets, are used because 
the series are shorter for the disordered systems but the 
triangular and rectangular shapes are still maintained. 
In Fig. [121 phase separation lines for three systems with 
the different degree of disorder characterized by the fol- 
lowing values of p — 0.5, p = 0.3 and p = 0.96 are shown. 
It can be clearly seen that Eq. (pij) describes the critical 
lines obtained by the PDA's very well around the homo- 
geneous point. For the p = 0.96 system, from fig = 1.5 
the critical points deviate from Eq. (pij) . Consistent with 
the results from the NPA's, the critical points obtained 
by PDA's deviate above the line given by Eq. ((3T|) . 

V. DISCUSSION & CONCLUSION 

To conclude, we have presented a detailed descrip- 
tion of the supercritical series expansions for the survival 
probability of the contact process in heterogeneous and 
disordered one-dimensional binary lattices. The hetero- 
geneous systems are modeled by lattices with repeating 
patterns of sites of two types A and B characterized by 
different recovery rates, fiA and hb, and the disordered 
systems are represented by lattices of similar sites ran- 
domly placed on the lattice sites with probabilities p and 
1 — p, for nodes A and J5, respectively. For the analy- 
sis of the two- variable series (in fiA and fis), we have 
presented a scheme based on NPA's in order to extract 
critical values and the critical exponent f3 and have also 
used PDA's to obtain estimates for the line of critical 
points. 

It has been demonstrated that (i) using symmetric and 
asymmetric linear transformations, it is possible to ex- 



tend the range of stable critical points in different regions 
of the rate-space (/iAj/is); (ii) keeping different orders in 
the Pade approximants with respect to the two variables, 
an extension of the line of critical points and an extended 
range for the critical exponents can be achieved; (iii) re- 
sults from NPA's and PDA's compare well, PDA's usu- 
ally widening the range of the critical lines a bit further 
in the (/i^, /iB)-plane; (iv) an approximate scheme for 
configurational averaging can be applied in disordered 
lattices. 

In general, the critical values can be reliably obtained 
by supercritical series expansions and they are in good 
agreement with the analytical approximation given by 
Eq. ((3T|) . For the critical exponents, the results are less 
conclusive due to larger errors, but certainly give some 
indication as to what happens to the universal critical be- 
havior when spatial heterogeneity and disorder is intro- 
duced. For all heterogeneous lattices, we see fluctuations 
away from the best known value for the CP from series 
expansions, /3 = 0.2769 [13], but they hardly ever exceed 
/3 = 0.28 over the range of reliable points. This suggests 
that the CP for heterogeneous lattices still belongs to the 
DP universality class. 

For the disordered systems, we see a qualitatively dif- 
ferent picture. In general, we find that the critical expo- 
nents monotonically increase away from the value at the 
homogeneous critical point. As we were only able to com- 
pute the exact configurationally averaged survival prob- 
ability up to order A'^ = 12 and an approximate one up 
to order N = 19, we do not have very high precision, but 
the above tendency is clearly visible in all the data. This 
picture of continuously changing exponents is consistent 
with what a number of other authors SBBHEil have 
found and with the Harris criterion [5|. 
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TABLE I: Series for the ultimate survival probability, Poo = 1 — 'l2n=i X]m=o CnmM^ "^Ms , for the heterogeneous lattice AB 
starting from a single occupied A-site up to order A'^ = 15 (continued in Tab. 
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TABLE II: Series for the ultimate survival probability, Poo = X^m=o '-"'"/^a ""Msi ^'^^ ^^e heterogeneous lattice AB 

starting from a single occupied A-site up to order A'^ = 24 (continued from Tab. |l]). 
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TABLE III: Series for the configurationally-averaged ultimate survival probability, (-Poo) = 1 — X)n=i ^m=oi'^nm) IJ-Ib, for 
the disordered lattice with impurity concentration p = 0.5 up to order N = 12. 
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